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Abstract. We explore a novel approach to the study of large-scale structure 
formation in which self-gravitating cold dark matter (CDM) is represented by a 
complex scalar field whose dynamics are governed by coupled Schrodinger and Poisson 
equations. We show that, in the quasi-linear regime, the Schrodinger equation can 
be reduced to the free-particle Schrodinger equation. We advocate using the free- 
particle Schrodinger equation as the basis of a new approximation method - the free- 
particle approximation - that is similar in spirit to the successful adhesion model. In 
this paper we test the free-particle approximation by appealing to a planar collapse 
scenario and find that our results are in excellent agreement with those of the 
Zeldovich approximation, provided care is taken when choosing a value for the effective 
Planck constant in the theory. We also discuss how extensions of the free-particle 
approximation are likely to require the inclusion of a time-dependent potential in the 
Schrodinger equation. Since the Schrodinger equation with a time-dependent potential 
is typically impossible to solve exactly, we investigate whether standard quantum- 
mechanical approximation techniques can be used, in a cosmological setting, to obtain 
useful solutions of the Schrodinger equation. In this paper we focus on one particular 
approximation method: time-dependent perturbation theory (TDPT). We elucidate 
the properties of perturbative solutions of the Schrodinger equation by considering a 
simple example: the gravitational evolution of a plane-symmetric density fluctuation. 
We use TDPT to calculate an approximate solution of the relevant Schrodinger 
equation and show that this perturbative solution can be used to successfully follow 
gravitational collapse beyond the linear regime, but there are several pitfalls to be 
avoided. 
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1. Introduction 

The local universe displays a rich hierarchical pattern of galaxy clustering that 
encompasses a vast range of length scales, culminating in rich clusters, super-clusters 
and filaments. However, the early universe was almost homogeneous with only slight 



Gravitational instability via the Schrddinger equation 



2 



temperature fluctuations seen in the cosmic microwave background radiation. Models of 
structure formation link these observations through the effect of gravity, relying on the 
fact that small initially over-dense regions accrete additional matter as the universe 
expands (a mechanism known as gravitational instability). The growth of density 
perturbations via gravitational instability is well understood in the linear regime, but 
the non-linear regime is much more complicated and generally not amenable to analytic 
solution. Numerical iV-body simulations have led the way towards an understanding 
of strongly developed clustering. Although such calculations have been priceless in 
establishing quantitative predictions of the large-scale structure expected to arise in a 
particular cosmology, it remains important to develop as full an analytical understanding 
as possible. After all, simulating a thing is not quite equivalent to understanding it. 

Analytical methods for studying the evolution of cosmological density perturbations 
fall into two broad classes (for a review, see pQ and references therein). First there are 
techniques based on applying perturbation theory (PT) to a hydrodynamical description 
of self-gravitating cold dark matter (CDM). These Eulerian approaches range from 
simple first-order (linear) PT (e.g. [2]) through to higher-order approaches of vastly 
increased complexity (e.g. [31 Hj). First-order Eulerian PT (the so-called linearized 
fluid approach) has been the mainstay of structure formation theory for many decades, 
by virtue of its simplicity and its robustness on large scales where density fluctuations 
are very much smaller than the mean density. When extrapolated to smaller scales, it 
provides a useful indicative measure of clustering strength, but can lead to absurdities 
(such as a negative matter density) if taken too far. Alternatively, there are Lagrangian 
approaches in which the trajectories of individual CDM particles are perturbed, rather 
than macroscopic fluid quantities (e.g. [3 El El E]). Even first-order Lagrangian PT, 
the celebrated Zeldovich approximation [9], is capable of following the gravitational 
collapse of density perturbations into the quasi-linear regime, beyond the breakdown 
of linear Eulerian PT. Comparisons with full iV-body calculations have shown that the 
Zeldovich approximation can successfully describe the quantitative morphology of the 
clustering pattern (e.g. \1Q\), as long as particle trajectories do not intersect. When 
particle trajectories cross (a phenomenon known as shell-crossing) the density field 
develops a formal singularity (known as a caustic). Since the Zeldovich approximation 
is purely kinematical, particles simply pass through the caustic and continue along their 
original trajectories. This leads to large-scale structure being rapidly 'washed out' in 
a manner not observed in iV-body simulations. A simple extension of the Zeldovich 
approximation that overcomes this problem is the adhesion approximation [UJ. In 
the adhesion approximation, particles move according to the Zeldovich approximation 
until trajectories intersect, but are assumed to 'stick' to each other when shell-crossing 
occurs. This sticking is achieved by including an artificial viscosity term in the equations 
of motion to mimic the action of the strong gravitational forces acting in the vicinity 
of a caustic. The adhesion model has proved very successful in explaining the pattern 
of large-scale structure observed in the universe; comparisons with iV-body simulations 
have shown impressive agreement far into the non-linear regime (e.g. [T2 | [T3 | [T H[T5l[T6] ). 
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In this paper we explore a radically different approach to the study of large- 
scale structure formation based on a wave-mechanical description of self-gravitating 
CDM. Following a suggestion by Widrow and Kaiser [17], we transform the usual 
hydrodynamical equations of motion into a non-linear Schrddinger equation coupled to 
a Poisson equation describing Newtonian gravity. The non-linear Schrddinger equation 
is similar to the familiar linear Schrddinger of quantum mechanics, except for the 
presence of an extra term known as the quantum pressure. The reasoning behind this 
transformation is not new - in essence it dates back to Madelung [18] - but it is quite 
neglected in cosmology. Coles [19] showed that the wave-mechanical approach provides a 
natural explanation of why the distribution function of density fluctuations, evolved from 
Gaussian initial conditions, should be close to the log-normal form that is observed. By 
studying the evolution of a one- dimensional sinusoidal density perturbation, Coles and 
Spencer [20] demonstrated that the wave-mechanical approach fares well in comparison 
with the Zeldovich approximation up until particle trajectories cross and that, when 
shell-crossing occurs, the wave-mechanical density field remains non-singular. Their 
work relied exclusively on numerical solutions of the Schrddinger equation; our objective 
here is to investigate how the Schrddinger-Poisson system can be solved by using 
approximation techniques instead. We show that, in an expanding universe, the 
Schrddinger equation can be reduced to the exactly solvable 'free-particle' Schrddinger 
equation by utilizing results from first-order Eulerian and Lagrangian PT. We advocate 
using the free-particle Schrddinger equation as the basis of a new approximation method 
- the free-particle approximation - that promises to be capable of evolving cosmological 
density perturbations into the quasi-linear regime. As we shall see, the free-particle 
approximation can essentially be thought of as an alternative to the adhesion model in 
which the viscosity term is replaced by the quantum pressure term. The performance of 
the new free-particle approximation is assessed, relative to the Zeldovich approximation, 
by considering a simple example of gravitational collapse. Particular attention is paid to 
elucidating the effect of the quantum pressure term since this has not been investigated 
previously. We also discuss how extensions of the free-particle approximation based on 
higher-order PT are likely to lead to a Schrddinger equation with a time-dependent 
potential that cannot be solved exactly. Consequently, we feel it is important to 
ascertain whether standard quantum-mechanical approximation techniques can be used, 
in a cosmological context, to obtain solutions of the Schrddinger equation with a time- 
dependent potential. In this work we concentrate on one particular approximation 
method: time-dependent perturbation theory (TDPT). Our goal is to investigate the 
properties of perturbative solutions of the Schrddinger equation in a cosmological setting; 
we perform our investigation by appealing to a simple gravitational collapse scenario. 
Although we restrict our attention to idealized examples in this work, the results we 
obtain will remain pertinent when the wave-mechanical method is applied to more 
general structure formation problems in the future. 

The outline of this paper is as follows: First, in order to keep the paper 
as self-contained as possible, we briefly outline the basic theory of cosmological 
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structure formation including the standard fluid approach, the Zeldovich approximation, 
the adhesion model, the wave-mechanical approach in general and the free-particle 
approximation in particular. Next, in section [3], we describe how the gravitational 
collapse of a plane-symmetric sinusoidal density perturbation can be modelled by using 
(i) the free-particle approximation and (ii) a perturbative solution of the Schrddinger 
equation with a time-dependent potential. We present our results in section 0] and 
conclude in section [5j 



2. Cosmological structure formation 

Observations suggest that there is approximately five times as much non-baryonic 
CDM in the universe as there is ordinary (luminous) matter. Consequently, large-scale 
structure formation is commonly studied by considering the gravitational amplification 
of fluctuations in the CDM distribution only. In standard treatments, collisionless CDM 
is assumed to be an ideal fluid with zero pressure (so-called dust) and the evolution 
equations for perturbations in the CDM fluid are obtained by linearizing the Einstein 
field equations about an expanding homogeneous and isotropic Friedmann-Robertson- 
Walker (FRW) cosmology (e.g. [2U E21 E3] ) . However, in this paper we will assume that 
the length scale of the perturbations is much smaller than the Hubble radius so that a 
Newtonian treatment is adequate (e.g. [2]). 

2.1. The fluid approach 

In an expanding universe, the Newtonian dynamical equations governing the evolution 
of fluctuations in a fluid of collisionless self-gravitating CDM can be written in the form: 

<9v 1 1 

_ + ffv + -(v.VJv + -V x * = 0, (1) 

£ + iv..[(l + *)v] = 0, (2) 
- 47rGa 2 pb,c<5 = 0, (3) 

where t is cosmological proper time and x = x(t) are comoving coordinates, related to 
physical coordinates r via r = ax. Here the scale factor a = a(t) has been normalized 
so that its value at the present epoch t is a = a(to) = 1. The Hubble parameter 
H = H(t) is defined by H = a/ a where a dot denotes a derivative with respect to t. 
The peculiar velocity field v = v(x, t) is given by v = ax and the potential $ = $(x, t) 
is the peculiar Newtonian gravitational potential. The density contrast 5 = 5(x,t) is 
5 = p/ph,c — 1 where p = p(x, t) is the CDM density field and pb,c = Pb,c(t) is the CDM 
density in the homogeneous FRW background. 

It is straightforward to show (e.g. [2]) from the fluid equations ([I]), (121) and ([31) that, 
to first-order in Eulerian PT, density perturbations in the CDM fluid grow according 
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to 5 = DSi where 8i is the density contrast at some initial time t\ and the linear growth 
factor D = D(t) is the growing mode solution of 

D + 2HD - AnGp htC D = 0, (4) 

normalized so that D\ = D(t{) = 1. The linear growth factor can then be used to write 
the fluid equations ([!]), (T5]) and ([3]) in the alternative form: 



f^ + (u-V x ) U+ ^u + V x e = 0, (5) 
f| + V x -[(l + 5)u] = 0, (6) 

v x e - ^ = o, (7) 

where the function / = f(D) is given by / = D/HD and Q c = Q C {D) is the familiar 
CDM density parameter: Q c = 8irGpb,c/3H 2 . The gravitational potential = 0(x, D) 
is related to the peculiar gravitational potential $ via = <&/a 2 D 2 . Similarly, the 
comoving velocity field u = u(x, D), defined by u = dx/dD, is obtained from the 
peculiar velocity v by a simple rescaling: u = w/aD. It is well known that linear 
Eulerian PT implies an irrotational velocity field: u = — V x 0, where the velocity 
potential = 0(x, D) is related to the gravitational potential via 



3a 



-0. (8) 



The linearized fluid equations only provide an accurate description of gravitational 
instability in the linear regime 5 < 1. A considerably more powerful method than the 
linearized fluid approach is the Zeldovich approximation [9] which is capable of following 
the evolution of cosmological density perturbations into the quasi-linear regime 5 ~ 1. 



2.2. The Zeldovich approximation 

The Zeldovich approximation is a Lagrangian approach (formally first-order Lagrangian 
PT) in which individual particle trajectories are considered: the comoving Eulerian 
coordinate x = x(q, D) of a particle initially located at a comoving Lagrangian 
coordinate q is given by 

x = q+(D-l)s, (9) 

where s = s(q) is a time-independent vector field. Requiring that mass is conserved 
immediately leads to the following expression for the CDM density field in the Zeldovich 
approximation: 

* = ^-i, do) 
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where S\ = 5i(q) is the initial density perturbation and J = J7"(q, D) is the determinant 
of the Jacobian of the mapping (Q between Lagrangian and Eulerian coordinates. The 
density field ( flOj) is known as the continuity density. 

The Zeldovich approximation can also be formulated in Eulerian space by noting 
that (E]) implies du/d.D = along a particle trajectory, provided particle trajectories do 
not intersect. Using the definition of the convective derivative d/dD = d/dD + u • V x , 
it follows that the Zeldovich approximation corresponds to 

— + (u-V x )u = 0, (11) 

in Eulerian space. Comparing (TlTI) with the Euler equation ([5]) it is apparent that the 
Zeldovich approximation also guarantees an irrotational velocity field u = — V x 0. We 
can then integrate (TTTT) to obtain the so-called Zeldovich- Bernoulli [21] equation: 

w ~ I = °- (12) 

which possesses an analytical solution = 0(x, D) of the form 

<P = <t>i-\(D-l)\s\\ (13) 

where fa = fa(q) is the initial velocity potential and x and q are related by the mapping 
([9]). As in the linearized fluid approach, the velocity potential <ft is related to the 
gravitational potential via (jSJ) and thus the Poisson equation (j7|) implies 

5 = DV x 0, (14) 

which is known as the dynamical density field. In general, the continuity and dynamical 
density fields are not equivalent. However, in one dimension, the density fields (II 01) and 
( fl4|) coincide and the Zeldovich approximation is an exact solution of the fluid equations 
as long as there is no shell-crossing (for a discussion, see p] and references therein). 

The Zeldovich approximation breaks down when particle trajectories cross since 
the mapping (Q from q to x is no longer unique. At a point where particle trajectories 
intersect, the velocity field becomes multi- valued, J = and the density field develops 
a formal singularity. The strong gravitational forces acting in the vicinity of the 
caustic should cause particles to be pulled towards it. However, since the Zeldovich 
approximation is kinematical, particles simply pass straight through the caustic and 
structure is rapidly 'smeared out' in an unphysical manner. The adhesion approximation 
[TT] was specifically designed to tackle this shortcoming of the Zeldovich approximation. 



2.3. The adhesion approximation 

In the adhesion approximation, particles follow Zeldovich trajectories until shell-crossing 
occurs. However, when particle trajectories cross, the particles are assumed to 'stick' to 
each other. As a result, the singularities predicted by the Zeldovich approximation 
are regularized and stable structures are formed, rather than being 'washed out'. 
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Mathematically, the adhesion model is obtained from the Zeldovich approximation by 
including an artificial viscosity term in ffTTj) so that: 



which is known as Burgers' equation |25j. The constant /i > has dimensions of L 2 
and can be thought of as a viscosity coefficient. A forcing term of the form V x ??, where 
r] = r](x, D) is a random potential, can be included in Burgers' equation to extend the 
adhesion model to the intergalactic medium (e.g. [26]). In the cosmologically relevant 
case of an irrotational velocity field u = — V x 0, (1151) can be integrated to obtain 



From a practical point of view, an attractive feature of Burgers' equation (or, 
equivalently, (|16l) ) is that it possesses an analytic solution. In the special case \i — > 
(known as the inviscid limit) a geometrical interpretation of the solution of Burgers' 
equation can be used to determine the 'skeleton' of the large-scale structure present 
at any given time (e.g. [TT| IT2| [15]). In this limit the structures formed in the 
adhesion model are infinitely thin and the adhesion approximation reduces exactly to 
the Zeldovich approximation outside of mass concentrations. For finite values of is, 
the viscosity term has an effect away from regions where particle trajectories cross and 
causes density perturbations to be suppressed on scales < /i 1//2 (e.g. [13J). The adhesion 
and Zeldovich approximations are then no longer identical outside of collapsing regions 
although, depending on the actual value of is, they become similar at a certain distance. 

The success of the adhesion model, relative to iV-body calculations, has led 
to adhesive gravitational clustering becoming an important concept in the study of 
cosmological structure formation. Although the original motivation for the adhesion 
approximation was purely phenomenological, it has recently been shown that, under 
certain simplifying assumptions, Burgers' equation can be naturally derived from the 
coarse-grained equations of motion [23[28j[29j[30l[3T]- In these modern interpretations of 
the adhesion model, the regularizing parameter \i emerges, not as a (constant) viscosity 
coefficient, but rather as a (density-dependent) gravitational multi-stream coefficient, 
arising from the self-gravitation of a multi-stream system. 

2-4- The wave-mechanical approach 

An alternative approach to the study of large-scale structure formation was suggested 
by Widrow and Kaiser [IT]. They proposed a wave-mechanical description of self- 
gravitating matter in which CDM is represented by a complex scalar field ip — i/)(pc, D) 
whose dynamics are governed by coupled Schrodinger and Poisson equations (see [32] 
for a relativistic extension of the original Newtonian theory). In this section we give a 
pedagogical derivation of the wave-mechanical formalism in an expanding universe since 
this has not been presented in the literature previously. 




(15) 




(16) 
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We begin by assuming an irrotational velocity field u = — V x 0. As discussed 
previously this is guaranteed to be the case in the linear and quasi-linear regimes and 
will remain so as long as there is no shell-crossing by Kelvin's circulation theorem. The 
Euler equation (jSJ) can then be integrated to obtain the Bernoulli equation 

w - \ ^ - v ^ < 17 > 

where the effective potential V = V(x, D) depends on both the gravitational potential 
G and velocity potential 4>: 

v = e-^. (is. 

The Bernoulli equation (fT7l) can be combined with the continuity equation ([6]) by 
performing a Madelung transformation [18] 



i>= (l + 5) 1/2 exp(^j, (19) 

where v is a real parameter with dimensions of L 2 . The wavefunction ip provides an 
elegant description of both density and velocity fields in a single complex function. 
Applying the Madelung transformation leads to the following coupled Schrodinger- 
Poisson system: 

i^=(-yV 2 +V + pV (20) 



V 



\) 3fic ( i \ 



W "' (N 2 -1)=0, (21) 



4V^ + VU. (23) 



2f 2 D 2 

where V = "P(x, D) is the so-called quantum pressure term, given by 
v 2 V 2 |^| 

p =y^r (22) 

If we neglect the quantum pressure term for the moment then we are left with the more 
familiar linear Schrddinger equation 

dip ( v 2 
dD 

Inserting the Madelung transformation into (1231) yields the usual continuity equation 
([6]) and a modified Bernoulli equation of the form 

f^|V x 0| 2 -V + P = O, (24) 

so the quantum pressure term appears in the fluid equations. This implies that we are 
free to drop the quantum pressure term from the non-linear Schrddinger equation ( 1201 ) 
and include it in the Bernoulli equation instead; this is the approach that we will adopt 
throughout this paper. 
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The parameter v appears in the Schrddinger equation in a manner analogous to the 
Planck constant in real quantum mechanics. However, the system we are considering 
here is entirely classical and thus v is treated as an adjustable parameter that controls 
the quantum pressure term V oc v 2 appearing in the Bernoulli equation (|24|) . The 
quantum pressure term acts as a regularizing term in the fluid equations, preventing the 
generation of multi-stream regions and singularities in the density field when particle 
trajectories cross. This was demonstrated by Coles and Spencer [20] who used the wave- 
mechanical approach to follow the gravitational collapse of a one- dimensional sinusoidal 
density perturbation through shell-crossing. The effect of the quantum pressure is thus 
qualitatively similar to that of the term fiV 2 ,^ in the adhesion model and v plays a 
similar role to the viscosity parameter /i, in the sense that they are both approximations 
to a general gravitational multi-stream coefficient arising in collisionless systems. The 
link between the wave-mechanical approach and the adhesion approximation will be 
further discussed below. In the semi- classical limit v — > we expect the effect of the 
quantum pressure term to be minimized and thus our wave-mechanical representation 
of self-gravitating CDM will approach the standard hydrodynamical description. 

The main deficiency of the wave-mechanical approach described here is that it is 
based on the fluid equations and thus implicitly assumes the existence of a single fluid 
velocity at each point. Consequently, the Schrddinger- Poisson system is incapable of 
properly describing the formation of multi-stream regions. However, this is not an 
intrinsic problem of the wave-mechanical formalism, but is instead a consequence of 
the simple Madelung form of the wavefunction we have used. Widrow and Kaiser [T7] 
showed that one can deploy more sophisticated representations of the wave function, 
such as the coherent-state formalism of Husimi [33], that allow for multi-streaming. 
The wave-mechanical approach can then be adapted to follow the full Vlasov evolution 
of the phase-space distribution function beyond the laminar flow regime [T71 1321 l3l] . 
These extensions are beyond the scope of this paper. 



2.4-1- The free-particle approximation The full Schrodinger-Poisson system cannot 
be solved analytically and so we must resort to either numerical or approximation 
methods; in this work we are concerned with the latter approach. The Schrodinger- 
Poisson system can be significantly simplified by using the fact that, to first-order in 
Eulerian and Lagrangian PT, the gravitational and velocity potentials are related via 
cf) = 2f 2 DQ/3Q c . In this case, it is evident from (f!8l) that the effective potential V 
is identically zero and the Schrddinger equation fl23|) reduces to the exactly solvable 
free-particle Schrddinger equation 

Accordingly, the Poisson equation f[2"TT) decouples from the Schrddinger equation and 
provides the following relationship between the amplitude and phase of the wavefunction: 

z/V^[arg(^)] + i(H 2 -l) =0. (26) 
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Inserting the Madelung transformation (1191) into (1251) and (J2H]) we obtain 

|^|Vx0| 2 + P = O, (27) 

along with the continuity equation (jH]) and the familiar relation (TH|) . The modified 
Bernoulli equation (|27|) looks like the adhesion approximation (fT6|) . except that the 
term has been replaced by the quantum pressure V . Since the quantum pressure 

term cannot be written in a form proportional to V^<P, we formulate an alternative 
to the adhesion model based on the free-particle Schrddinger equation (1251) and the 
Poisson-like equation (|26|) . We call this new approximation scheme the free-particle 



approximation. In Appendix A we show that, in the semi-classical limit v — > 0, the free- 
particle approximation reduces to the Zeldovich approximation prior to shell-crossing. 
The quantum pressure term will then only become important in regions where particle 
trajectories intersect. 

3. Examples of gravitational collapse using the Schrodinger equation 

3.1. A simple test of the free-particle approximation 

The free-particle approximation should be capable of following the evolution of density 
fluctuations into the quasi-linear regime. In order to verify this hypothesis, we test the 
free-particle approximation by appealing to a simple example: the gravitational collapse 
of a plane-symmetric density perturbation. The assumption of plane-symmetry implies 
that we are effectively considering a one-dimensional problem. This is advantageous 
since, in one dimension, the Zeldovich approximation provides an exact solution of the 
equations of motion until shell-crossing occurs. The performance of the new free-particle 
approximation can then be assessed relative to the exact Zeldovich solution. We will 
consider the case where the parameter v is finite (rather than v — *• 0) and investigate the 
effect of the quantum pressure term as v is varied. This is an important issue to address 
since the quantum pressure term is an integral part of the free-particle approximation 
(and, indeed, the wave-mechanical approach in general). For illustrative purposes we 
assume the plane-symmetric initial density perturbation b~\ = 5i(x) is sinusoidal: 

4 = ^008 ^), (28) 

where d is the comoving period of the perturbation and < |5 a | < 1 to ensure that the 
initial CDM density field is everywhere non-negative. The corresponding initial velocity 
potential perturbation is found from (TSHJ): 

* = -(£)'*■ (29) 

The one-dimensional free-particle Schrodinger equation 
. dip v 2 d 2 ^ 

w w = -y^ (30) 
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must then be solved subject to the initial condition 

V* = (l + ft)^expf^\ (31) 

with Si and 0; given by (128]) and (1291) . respectively. In order to simplify calculations 
we consider a large cubic volume of comoving side length L equipped with periodic 
boundary conditions at each face. This is a construction commonly used in the study of 
cosmological structure formation; the limit L — > oo can always be taken as a final step. 
Furthermore, we divide the cubic volume into cells of side length d (i.e. we set L = Nd, 
N > an integer) since the initial density perturbation (1281) is periodic with comoving 
period d. 



3.1.1. Solution of the free-particle Schrddinger equation The solution of the free- 
particle Schrddinger equation in a static background cosmology is discussed in detail 
in Appendix B.l Rather than repeating the calculation for the case of an expanding 



background, it suffices to note that we can simply replace proper time t by the linear 
growth factor D in the solution obtained in Appendix B.l The solution of the free- 
particle Schrddinger equation fl30l) is then 



a n exp 



-i(D - 1)ES 



(o) 



(32) 



where En^ 
eigenfunctions 



v 2 k 2 j2 and k n 
f = 4 0) (x) 



= 2nir/d is a comoving wavenumber. The (orthonormal) 
are of the form 



L 3 / 2 



exp (ik n x) 



(33) 



and the expansion coefficients determined from 



NL 2 



4>n\x)lpi(x)dx, 



(34) 



where the over-line denotes complex conjugation and the initial wavefunction t/>; is given 
by fl3TT) . It follows from fl33|) that the expansion coefficients a n can be calculated simply 
by taking the Fourier transform of the initial wavefunction. 



3.2. Beyond the free-particle approximation 

The free-particle formalism provides an approximation to the full Schrodinger-Poisson 
system that is exactly solvable. The free-particle approximation is motivated by the 
fact that the effective potential V appearing in the Schrodinger-Poisson system is zero 
to first-order in both Eulerian and Lagrangian PT. Since the free-particle approximation 
is based on linear PT, a natural extension of the free-particle approximation could be 
constructed simply by including higher-order terms in PT. However, including higher- 
order terms in PT implies V ^ in general (e.g. [31 [8j), leading to a time-dependent 
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external potential in the Schrddinger equation (I23p . The Schrddinger equation with 
a time-dependent potential is typically impossible to solve analytically, so either a 
numerical or an approximate solution must be sought. 

In light of the above comments, we feel it is important to investigate whether 
quantum- mechanical approximation methods can be used, in a cosmological context, to 
obtain useful solutions of the Schrddinger equation with a time-dependent potential. 
We will focus on time-dependent perturbation theory (TDPT) since this is one of 
the most widely used approximation schemes in quantum mechanics. Our aim is to 
elucidate the properties of perturbative solutions of the Schrddinger equation with a 
time-dependent potential in a cosmological setting. In particular, we wish to assess 
the effect of varying the parameter v where we again assume that v is finite. It is not 
our intention here to explicitly formulate a higher-order extension of the free-particle 
approximation; we leave this for future work. Instead we investigate perturbative 
solutions of the Schrddinger equation by considering a simple idealized scenario. We 
begin by rewriting the Schrddinger-Poisson system in a static background cosmology; 
the reason for this will become apparent. In a static universe the scale factor is constant 
a = oq = 1 (i.e. H = 0) and the Friedmann equations imply that the universe is closed 
with a non-zero cosmological constant A = 47rGpb,c- The Schrddinger-Poisson system 
can then be written in the form 



■4H-t v H^ (35) 

V^-4vrGX c (H 2 -l) =0, (36) 

where comoving coordinates x now coincide with physical coordinates r. The potential 
$ is the peculiar gravitational potential and 

V>=(l + 5) 1/2 exp(^, (37) 

where the peculiar velocity v is the gradient of the velocity potential (p = y(x, t): 
v = — V x y3. Inserting the Madelung transformation (1371) into (I35I) and (|36|) we obtain 

^-I|V x y?| 2 -$ + P = 0, (38) 

as well as the continuity equation (T5]) and the Poisson equation Q in a static universe. 
The gradient of the modified Bernoulli equation ( 1381) gives the Euler equation (CQ) in a 
static background, but with an extra term corresponding to the gradient of the quantum 
pressure V . It follows from the fluid equations (CD), (TSJ) and ([3]) in a static background 
that, to first-order in Eulerian PT, the peculiar gravitational potential evolves according 
to $ = -D^i, where D = exp \(t — t;)/r] and r = 1/A 1 / 2 is the characteristic time-scale 
for the collapse of a density fluctuation. To first-order in Lagrangian PT (i.e. the 
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Zeldovich approximation) it is straightforward to show from (jSJ) and (I13p that, in a 
static universe, the gravitational potential $ = $(x, t) is given by 



$ = D 



$i-^(£>-l)N 



(39) 



where = $i(q) is the initial gravitational potential, Vj = Vi(q) is the initial 
peculiar velocity field and the coordinates x and q are related by the mapping ([9]) 
with s = tv;. Therefore, in a static background, we have a non-zero time-dependent 
potential appearing in the Schrddinger equation (1351) to first- order. This is in contrast 
to the expanding background case where we must include higher-order terms in PT 
to obtain V ^ in the Schrddinger equation (|23|) . Since our goal is to investigate 
perturbative solutions of the Schrddinger equation with a time-dependent potential, it 
will clearly be simpler to work in a static background as we can then use a first-order 
(rather than a cumbersome higher-order) expression for the potential. This is why we 
have assumed a static background cosmology. 

To proceed with our investigation, we now apply the wave-mechanical approach to 
a specific example of gravitational collapse. As before, we consider the evolution of a 
plane-symmetric sinusoidal density perturbation b~\ of the form (}28l) . The corresponding 
initial peculiar gravitational potential is then found from the Poisson equation (136|) : 

$ i = "f^-V^ (4°) 



2itt t 

and the initial peculiar velocity field is given by 

9 — Y (<£"*?), (41) 

which follows from (j40l) by noting that, in a static universe, the linear (Eulerian and 
Lagrangian) PT relation ([HD implies the gravitational potential $ and velocity potential 
<p satisfy ip = r$. The full one-dimensional Schrddinger equation we are aiming to solve 
is then 

dip f v 2 d 2 
~dt ~ \^2dx 2 

where we use the first-order Lagrangian PT expression (13"9~1) for the gravitational 
potential $ = $(x, t) since this is exact in one dimension (up until shell-crossing occurs). 
The initial gravitational potential $; = and peculiar velocity v\ = v\(q) in (1391) 

are given by (HOI) and (pfll) . respectively. The Schrddinger equation (H2|) must be solved 
subject to the initial condition 

^(l + aO^expf^Y (43) 



^ = ( ~~ + * ) (42) 



where 5{ is given by (1281) and <p\ = r$i. We again restrict our attention to a large cubic 
volume of side length L equipped with periodic boundary conditions at each face. As 
before, the periodicity of the initial density perturbation suggests it will be convenient 
to partition the cubic volume into cubic cells of side length d by writing L = Nd, N > 
an integer. 
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3.2.1. Perturbative solution of the Schrddinger equation with a time- dependent potential 
The Schrddinger equation (|42p cannot be solved analytically; we have used TDPT to 
determine a second-order approximate solution instead. The details of the calculation 



are presented in Appendix B.2 We summarize the main results here for reference. The 
perturbative solution was found to be 



3=0 



(44) 



where the zeroth-order term ip^ is simply the solution of the free-particle Schrddinger 
equation ($ = 0): 

~ -i(t-fi)£&° r 



■0 



(0) 



a n exp 



(45) 



with En = u 2 k 2 l /2 and k n = 2nir/d. The eigenfunctions (ffi are of the form ( 1331) and 
the expansion coefficients a n are found from (1341) with the initial wavefunction ip\ now 
given by ( 1431) . The first-order term i s 



(i) 



4 



(46) 



n m 



where 7 = vt jd 2 is a dimensionless parameter, Z~ TO>n = X mn {t) is defined by 



1 /■* 

Im,n = - / dt exp 



-i(t - t')E, 



(o) 

m 



$ m , n (t')exp 



-i(f - fi)£, 



(o) 

n 



(47) 



and the matrix elements $ m>n = are given by 



J p (2pa) + - ^ J p+S (2pa) 



s=±2 



+ (y - x ) E Vs(2p«) 



(48) 



=±3 



where p = m — n is an integer, a = a(t) is defined by a = 5 a (D — l)/2 and the functions 
Jz are Bessel functions of the first kind. The second-order term ijj^ is given by 



1> 



(2) 



(0) 



(49) 



m I 



where /Q )mn = £j, m , n (i) is defined by 
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K 



Lm,n 



T- 



dt' / dt"exp 



x$ m>ri (t") exp 



-i(t-f)El 



(o) 



V 



$ z , m (t')exp 



7/ 



-i(t" - 



(0) 



(50) 



In order to complete the perturbation expansion (|44p of the wavefunction ip, the integrals 
f|47|) and fl50|) must be evaluated numerically. 



4. Results and discussion 

We now analyze the gravitational collapse of the plane-symmetric sinusoidal density 
perturbation ( 1281) in (i) an expanding background, using the free-particle approximation, 
and (ii) a static background, using a perturbative solution of the Schrddinger equation 
with a time-dependent potential. In both cases the CDM density field is obtained from 
the amplitude of the wavefunction via 5 = \ip\ 2 — 1 and compared with the linear growth 
law 5 = D5\ and the exact Zeldovich solution 

(1 + *0 



S 



(51) 



[1 — (D — 1)4] ' 

which follows from (TTUT) and (|T4|) . The evolution of the density fluctuation is followed up 
until particle trajectories cross which, from (15T|) . occurs when the linear growth factor 
D = D sc = 1 + l/<5 a . Hereafter we will assume, for illustrative purposes, that the 
amplitude of the initial density perturbation is 5 a = 0.01, so that D sc = 101. 



Planar collapse using the free-particle approximation 

The (exact) solution of the one- dimensional free-particle Schrddinger equation (1301 is 
given by (|32|) . Upon introducing the dimensionless comoving coordinate x = x/d, 
we find that the free-particle solution depends solely on the dimensionless parameter 
T = v/d 2 , which we assume to be finite (as opposed to V — > 0). Recall that, in the 
semi-classical limit r — > 0, the free-particle approximation coincides with the Zeldovich 
approximation before shell-crossing occurs (see Appendix A ). In one dimension this 
implies that the free-particle approximation will be exact in the limit r — > (provided 
there is no shell-crossing). This suggests that our numerical implementation of the free- 
particle approximation with finite T can be optimized by using the smallest possible 
value of T. In practice there is a lower bound T c on the value of T arising from the fact 
that we test the free-particle approximation on a discrete grid. Sampling at the Nyquist 
rate (the minimum possible sampling rate) requires that the phase change between 
two neighbouring grid points must be less than or equal to n radians. The phase of the 
initial wavefunction is proportional to 1/r and so, as V approaches zero, the phase varies 
increasingly rapidly and the change in phase between two neighbouring grid points can 
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exceed it radians. The phase is then insufficiently sampled and aliasing effects cause 
our free-particle method to break down. To avoid this problem we must choose T > r c ; 
we expect the performance of the free-particle approximation to be optimal for T = T c . 
It is straightforward to show from (l3Tj) that r c = <5 a /27r 2 iVg, where N g is the number 
of grid points. We use N g = 512 so that T c = 1 x 1CT 6 . We now investigate how the 
free-particle approximation, and in particular the quantum pressure term, behaves as V 
is decreased from a large value towards T = r c . 

We begin by assuming T = 1. The plots in the left-hand column of figured] show 
the evolution of the free-particle density field 5 = \ip\ 2 — 1 in this case. The density 
field is only plotted in the interval < x < d since, by construction, the density field is 
periodic with comoving period d. We can immediately see that the free-particle density 
field simply oscillates about the mean value (5) = and there is no net growth of the 
initial density perturbation. We find that this oscillatory behaviour also persists for any 
T > 1. This effect is caused by the quantum pressure term. Recall that the quantum 
pressure V enters the free-particle formalism via the modified Bernoulli equation ( 1271) . 
Since V oc v 2 oc T 2 , the parameter T controls the size of the quantum pressure term 
relative to the convective term C = C(x, D), defined by C = — |V x 0| 2 /2. The plots 
in the right-hand column of figure [1] compare the magnitudes of the quantum pressure 
and convective terms as a function of the linear growth factor. The quantum pressure 
term is the dominant term appearing in ( [271) at all times (except when it momentarily 
becomes zero as the density field passes through 5 = 0). Therefore, the growth of the 
density perturbation appears to be inhibited by the large quantum pressure term. This 
is somewhat reminiscent of the adhesion approximation where the growth of fluctuations 
is suppressed on scales < yU 1 / 2 , leading to the large-scale structure distribution appearing 
'washed out' for large [i [13]. However, even if the viscosity term is large, the density 
field in the adhesion model does not display oscillatory behaviour; this is a peculiarity 
of the quantum pressure term. It is interesting to note that gas pressure in a baryonic 
fluid causes a qualitatively similar effect (e.g. [2]). If we were considering a baryonic 
fluid (rather than a pressureless fluid of CDM) then it would be necessary to account for 
the effects of gas pressure by including a classical pressure term in the Euler equation. 
Applying first-order Eulerian PT to the fluid equations then yields a characteristic length 
scale (the Jean's length)] any Fourier modes of the density field with a wavelength smaller 
than the Jean's length undergo damped oscillation rather than growth. However, we 
must be careful not to directly compare the effects of the quantum pressure term to 
those of a classical pressure term since the origin and form of the two terms is very 
different. 

The quantum pressure V oc T 2 should become less significant as T is decreased. 
The left-hand column of figure [2] displays the free-particle density field as a function 
of the linear growth factor when T = 1.4 x 10~ 3 . Figure EJ^a) and figure E](b) show 
that the free-particle density field evolves as expected up to D = D sc /2, in the sense 
that matter flows towards initially over-dense regions. However, the growth rate of the 
over- densities is considerably less than that predicted by the exact Zeldovich solution. 
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Figure 1. Evolution of a plane-symmetric sinusoidal density perturbation in an 
expanding CDM-dominated universe. The amplitude of the initial density fluctuation 
is <5 a = 0.01. The left-hand plots show the density fields obtained from the free- 
particle approximation (FP), the linearized fluid approach (LEA) and the Zeldovich 
approximation (ZA), at three different values of the linear growth factor. The 
parameter T = 1 in the free-particle approximation. The right-hand plots show the 
corresponding evolution of the magnitudes of the convective (C) and quantum pressure 
(P) terms, in units of d 2 . 



For D > D sc /2, the over-densities in the free-particle density field cease to grow and 
actually begin to decay, leading to a large discrepancy between the free-particle and 
Zeldovich approximation density fields at D ~ D sc ; see figure 12(c). This behaviour is 
counter-intuitive; over-densities should continue to grow up until shell-crossing occurs 
as in first-order Eulerian and Lagrangian PT. The quantum pressure term V again 
provides the explanation. The right-hand column of figure [2] shows the evolution of \C\ 
and \V\. We can see that the convective term initially dominates the quantum pressure 
term in the Bernoulli equation (|27|) . This remains the case at early times, corresponding 
to when the free-particle approximation agrees well with the linearized fluid approach 
and the Zeldovich approximation. However, the quantum pressure quickly grows until, 
at D « D sc /2, it is much larger than the convective term. The large quantum pressure 
term suppresses the gravitational collapse of the density perturbation, causing the over- 
densities to cease growing and subsequently decay. 
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Figure 2. Evolution of a plane-symmetric sinusoidal density perturbation in an 
expanding CDM-dominated universe. The layout of the plots is as in figure[IJ the only 
difference is that the parameter r = 1.4 x 10~ 3 in the free-particle approximation. 



The effect of the quantum pressure term is indeed less pronounced in the T = 
1.4 x 1CT 3 case than in the T = 1 case. However, at late times it causes matter to 
flow away from over-dense regions in an unrealistic manner. We find that, when F 
is reduced to T « 1 x 10~ 3 , the quantum pressure never becomes large enough to 
cause the decay of over-densities and the initial density perturbation grows up until 
shell-crossing. However, the collapse of the fluctuation is still highly suppressed by the 
quantum pressure term in high-density regions and the free-particle approximation thus 
provides a poor match to the Zeldovich approximation in such regions. As T — > r c , the 
suppression effect of the quantum pressure term diminishes and the agreement between 
the free-particle and Zeldovich approximations improves. The left-hand column of figure 
[3] shows the free-particle evolution of the initial density perturbation in the optimal 
case T = T c = 1 x 10~ 6 . At all values of the linear growth factor the free-particle 
approximation now provides an excellent match to the exact Zeldovich solution, even 
in over-dense regions. However, unlike the Zeldovich approximation, the free-particle 
approximation leads to a density field that remains well behaved at shell-crossing. The 
right-hand column of figure [3] shows the corresponding evolution of \C\ and \V\. It 
is immediately apparent that the convective term completely dominates the quantum 
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pressure term at all values of the linear growth factor. Although the quantum pressure 
term is small, it is non-zero (because we are using a finite T) and thus will have some 
effect on the collapse of the density perturbation. However, this effect is negligibly small 
for T = 1 x 1(T 6 . 
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Figure 3. Evolution of a plane-symmetric sinusoidal density perturbation in an 
expanding CDM-dominated universe. The layout of the plots is as in figure[ll the only 
difference is that the parameter f = f c = 1 x 10~ 6 in the free-particle approximation. 



4-2. Planar collapse using the Schrddinger equation with a time- dependent potential 

In a static background the linear growth factor is given by D = exp [(t — t\)/r\. Upon 
introducing a dimensionless time coordinate t — (t — t\)/r we find that, for 5 a = 0.01, 
shell-crossing occurs when i = i sc = In (D sc ) m 4.6. We now use the perturbative solution 
( |44|) of the one- dimensional Schrddinger equation with a time-dependent potential 
( 1421) to follow the gravitational collapse of the initial density perturbation ( 1281) up 
to i sc . Observe that the perturbation expansion of the wavefunction depends on one 
dimensionless parameter: 7 = ur/d 2 , which we assume to be finite. As in the previous 
section, phase-aliasing effects impose a lower bound 7 C on the value of 7 in practice. It 
follows from (|43j) that again 7 C = S Su /2n 2 N g = 1 x 10~ 6 . Based on the results of our 
test of the free-particle approximation, it is reasonable to suppose that the performance 
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of the wave-mechanical approach will be optimized for 7 = 7 C - Since the gravitational 
potential fl39l) appearing in fj42|) is exact in one-dimension, we would expect the wave- 
mechanical density field to agree well with that of the exact Zeldovich solution in this 
optimal case. Unfortunately, as we shall see, the situation is more complex when dealing 
with perturbative solutions of the Schrddinger equation. We now examine the behaviour 
of the TDPT solution of the Schrddinger equation as 7 is reduced from a large value 
towards 7 = 7 C - 

As a starting point we choose 7 = 1. In this case we find that the wave- mechanical 
density field 5 = IV'I 2- 1 oscillates rapidly about (5) = and there is no overall growth of 
the initial density fluctuation (cf. figure 1). In fact, the wave-mechanical approach leads 
to an oscillatory density field for any 7 > 1. As in the previous section, this is due to the 
fact that the quantum pressure term is the dominant term in the Bernoulli equation ( |38l ) 
for all times up to shell-crossing. The quantum pressure V oc v 2 oc 7 2 and thus should 
become less influential as 7 is decreased. The evolution of the wave-mechanical density 
field for 7 = 0.1 is shown in the left-hand column of figure HI By comparing figure H^a) 
with figure IU(b) we see that the initial over-densities begin to decay rather than grow; 
this is unphysical. However, at t ~ i sc /3, the decay slows and halts and matter begins 
to move towards over-dense regions in the usual fashion up to t ~ i sc . Since over-dense 
regions only begin to collapse after i ~ i sc /3, the final wave-mechanical density field 
is considerably more homogeneous than the density fields obtained from the linearized 
fluid approach and the Zeldovich approximation. The reason for this behaviour becomes 
apparent if we examine the relative sizes of the convective, gravitational potential and 
quantum pressure terms in (|38|) . The plots in the right-hand column of figure H] compare 
the time- variation of \C\, |$| and \V\. At i = the quantum pressure is the dominant 
term and subsequently causes the initial over-densities to decay. However, in a static 
background, the gravitational potential (139]) grows rapidly with time and, at i w t sc /3, 
dominates over the quantum pressure and convective terms. This corresponds to the 
time when the over-densities cease to decay and the gravitational collapse of the density 
fluctuation begins to proceed as expected. The gravitational potential then remains the 
dominant term up until shell-crossing occurs. 

We find that, as the value of 7 is reduced further, the effect of the quantum pressure 
term indeed lessens. However, for 7 < 0.02, the wave-mechanical approach exhibits 
rather unusual behaviour. For example, consider 7 = 0.015 (note 7 ^> 7 C ); figure 
[5] illustrates the behaviour of the wave-mechanical approximation in this case. It is 
evident from figure [5](a) and figure [S^b) that the wave-mechanical approximation now 
provides a close match to the linearized fluid approach and the exact Zeldovich solution 
up to i Rj 2t sc /3. However, at times close to shell-crossing, an over-density appears at 
x = d/2 (i.e. where the gravitational potential is at a maximum) in the wave-mechanical 
density field. This is clearly not realistic. To explain this, recall that the second-order 
perturbative solution of the Schrddinger equation fH2|) is 
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Figure 4. Evolution of a plane-symmetric sinusoidal density perturbation in a 
static CDM-dominated universe. The amplitude of the initial density perturbation 
is 6 a = 0.01. The left-hand plots show the density fields obtained from the wave- 
mechanical approach with a time-dependent potential (WM), the linearized fluid 
approach (LFA) and the Zeldovich approximation (ZA). The parameter 7 = 0.1 in 
the wave-mechanical approximation. The right-hand plots show the corresponding 
evolution of the magnitudes of the convective (C), gravitational potential (<&) and 
quantum pressure (P) terms, in units of d 2 /t 2 . 



2 

V = J> W , (52) 
3=0 

where ip(°\ ip^ and ip^ are given by (j45p . (j46p and (j49p . respectively. In order for the 
perturbation expansion fl52|) to be valid we require the first and second-order terms to 
be much smaller than the zeroth-order term. However, ip^ 1 ' oc I/7 and i/j^ oc I/7 2 so, 
for small values of 7, we anticipate that these terms will be large compared to ip(°\ To 
check this we define a ratio 

- 

with j = 1, 2. We then require xi,o *C 1 and x 2 ,o 1 for our perturbative solution 
to hold. The ratios xi,o an d X2,o are shown as functions of time in the right-hand 



Gravitational instability via the Schrddinger equation 



22 




0.2 0.4 0.6 0.8 1 0.2 0.4 0.6 0.8 1 

x/d x/d 



Figure 5. Evolution of a plane-symmetric sinusoidal density perturbation in a static 
CDM-dominated universe. The left-hand plots are as in figure 01 except that the 
parameter 7 = 0.015 in the wave-mechanical approximation. The right-hand plots 
show the corresponding evolution of the ratios X1.0 and \2,a defined in the text. 

column of figure At early times xi,o an d X2,o are both small as they should be. 
However, they both grow with time until, at t ~ i sc , the first and second-order terms 
in the perturbation expansion of the wavefunction are comparable to the zeroth-order 
term. A general feature is that the first-order term initially grows at a faster rate than 
the second-order term but, as the time of shell-crossing approaches, the second-order 
term begins to grow very rapidly and is soon comparable to the first-order term. The 
false over-density at x = d/2 then develops, suggesting that this is due to the onset of 
non-perturbative behaviour in the system. We find that, as 7 — > 7 C , the perturbative 
solution of the Schrddinger equation breaks down at progressively earlier times. 

The development of non-perturbative behaviour is somewhat reminiscent of the 
onset of turbulence in a classical fluid. It is interesting to speculate that small values of 
7 could be linked to some form of turbulent behaviour in our hydrodynamical description 
of quantum mechanics, especially when we realise that v has dimensions of viscosity and 
that the dimensionless parameter I/7 is similar to the Reynold's number appearing in 
classical fluid mechanics. Turbulence in classical fluids is associated with the formation 
of vortices which obviously does not occur in one-dimensional systems. Nevertheless the 
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analogy is compelling. A possible explanation is that very small values of 7 involve large 
angular frequencies in the phase of the wavefunction (which is proportional to I/7). To 
regain classical behaviour in the limit 7 — >■ the waves from which tp is constructed must 
cancel exactly. Since the limit 7^0 involves waves of infinitely small wavelength, any 
finite perturbative calculation will find convergence impossible. A related but converse 
phenomenon arises in Eulerian fluids (i.e. fluids without intrinsic viscosity) when they 
are described in terms of a truncated spectrum [35] . 

It is clear from the above discussion that there is a conflict to be faced when 
choosing an appropriate value of 7. On the one hand, we wish to set 7 = 7 C to ensure 
that the quantum pressure term is minimized. On the other hand, the accuracy of our 
perturbative solution of the Schrddinger equation improves as 7 is increased. The best 
one can do in this situation is to compromise by selecting the smallest possible value 
of 7 for which the perturbative solution of the Schrddinger equation remains physically 
reasonable over the time-scales of interest. We denote this 'optimal' value by j Q . In 
the situation under consideration here, we find that 7 Q = 0.02; the left-hand column of 
figure [6] shows the behaviour of the wave-mechanical density field in this case. We can 
see from figure E](a) and figure Mjo) that the wave- mechanical approach is again in good 
agreement with the linearized fluid approach and the exact Zeldovich solution for times 
up to t ps 2i sc /3. Figure [6](c) shows that, at i ps t sc , the wave-mechanical density field 
is well-behaved (in the sense that there is no unphysical over-density at x — d/2) and 
agrees extremely well with the Zeldovich approximation in the vicinity of x = d/2. The 
wave-mechanical approach also leads to the formation of over-densities with 5 ~ 2.2 
(i.e. in the quasi-linear regime). However, there is still a large discrepancy between the 
wave-mechanical and Zeldovich approximation density fields in high-density regions. 
This is partly attributable to the fact that, although dominated by the gravitational 
potential term, the quantum pressure term is still non-negligible for 7 = 0.02 and thus 
hinders the collapse of the density fluctuation. In addition, as shown by the plots in 
the right-hand column of figure [6j the ratios Xi,o an d X2,o grow with time, eventually 
becoming comparable to unity at i ~ £ sc . Consequently, our perturbative solution (1441) 
of the Schrddinger equation is on the verge of failing completely at times close to shell- 
crossing. 

5. Conclusion 

In this paper we have proposed a new approximation method suitable for evolving 
cosmological density perturbations into the quasi-linear regime. This new method 
is based upon the free-particle Schrddinger equation and is thus called the free- 
particle approximation. We have seen that the free-particle approximation provides 
an alternative to the well-known adhesion model in which the viscosity term /iV^ is 
replaced by the quantum pressure term V . The quantum pressure acts as a regularizing 
term that prevents the generation of density singularities when particle trajectories cross. 
Therefore, in contrast to the Zeldovich approximation, the free-particle method does not 
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fail catastrophically when shell-crossing occurs. The quantum pressure term is controlled 
by a free parameter v which plays a similar role to the viscosity parameter //. We have 
shown that, in the semi-classical limit v — * 0, the free-particle approximation reduces to 
the Zeldovich approximation, provided there is no shell-crossing. The quantum pressure 
term will then only have an effect in multi-stream regions. 

We have performed a simple test of the free-particle approximation by considering 
the gravitational collapse of a plane-symmetric sinusoidal density perturbation. The 
behaviour of the free-particle approximation has been examined for several different 
finite values of the parameter v (or, equivalently, T in the notation of section 14. ip in 
an attempt to elucidate the role of the quantum pressure term. We have found that, 
contrary to the v — > case, the quantum pressure affects the free-particle approximation 
before shell-crossing occurs if v is finite. In particular, the quantum pressure acts to 
suppress the gravitational collapse of density fluctuations; the degree of suppression 
depends on the actual value of v. For large values of u, the quantum pressure term 
completely dominates over the convective term in the modified Bernoulli equation ( 1271) 
and the free-particle density field then simply oscillates about (S) = 0. As the value of 
v is decreased the effect of the quantum pressure term diminishes. When the value of v 
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is set to the smallest possible value u c allowed by our numerical implementation of the 
free-particle method, we find that the free-particle approximation provides an excellent 
match to the exact Zeldovich solution of the equations of motion, up until shell-crossing 
occurs. In this case the quantum pressure term is still non-zero, but is negligibly small 
in comparison with the convective term and thus has little effect before shell-crossing. 

The success of the free-particle approximation in the simple example presented 
here suggests that this method promises to be a useful analytical tool for modelling 
the formation of large-scale structure in the universe. In a subsequent paper we will 
exhaustively test the free-particle approximation in a more cosmologically relevant 
scenario by appealing to a full iV-body simulation |36j . 

The free-particle approximation is motivated by the fact that, to first-order in 
Eulerian and Lagrangian PT, the effective potential V appearing in the full Schrddinger- 
Poisson system is zero. In future work it may prove desirable to go beyond the free- 
particle approximation by including higher-order terms in PT. We anticipate that any 
such extension of the free-particle approximation will, in general, lead to a Schrddinger 
equation with a time-dependent potential that cannot be solved exactly. Instead, either 
numerical or approximation methods must be deployed. In this work we have examined 
whether TDPT is a suitable technique for constructing approximate solutions of the 
Schrddinger equation in a cosmological setting. We have performed our investigation 
by appealing to a simple model where the universe is assumed to be static. In this case 
the gravitational potential $ appearing in the Schrodinger-Poisson system is non-zero 
to first-order in both Eulerian and Lagrangian PT. We use first-order Lagrangian PT 
(which is exact in one-dimension) to calculate $ in the Schrddinger equation for the 
case of a plane-symmetric sinusoidal density perturbation. A perturbative solution of 
the resulting Schrddinger equation was then found by applying second-order TDPT. To 
elucidate the properties of this perturbative solution we have studied the behaviour of 
the wave-mechanical approximation as the (finite) parameter v (or 7 in the notation 
of section 14. 2[) is varied. For large values of v the quantum pressure term V is the 
dominant term in the modified Bernoulli equation (1381) and again acts to inhibit the 
gravitational collapse of the density fluctuation. On the other hand, we have seen 
that the perturbation expansion of the wavefunction breaks down for small values of u, 
leading to the formation of an unphysical over-density where the gravitational potential 
is at a maximum. Therefore, in order to optimize the performance of our perturbative 
wave-mechanical approach, we must be satisfied with choosing the smallest possible 
value of v for which the TDPT solution of the Schrddinger equation remains physically 
sensible over the time-scales of interest. In this optimal case we find that the wave- 
mechanical approach is capable of evolving the density perturbation into the quasi-linear 
regime, but provides a poor match to the exact solution in over-dense regions. 

In summary, we have shown that TDPT is a viable means of calculating 
approximate solutions of the Schrddinger equation in a cosmological context, provided 
extreme care is taken in choosing v. Although we have only considered a simple example, 
we believe that the issues we have highlighted regarding perturbative solutions of the 
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Schrddinger equation will remain pertinent in more general situations. Thus our work 
promises to be useful in the future when constructing more sophisticated alternatives 
to the free-particle approximation. 
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Appendix A. Path integral solution of the free-particle Schrodinger 
equation in the semi-classical limit 

The solution ip = ipfa, D) of the free-particle Schrodinger equation 
• dip v 2 2 

can be written in the form 

V = / C/(x, J D;q,l)Vi(q)d 3 q, 



(A.l) 



(A.2) 



where ip\ = ipi(q) is some initial wavefunction and U = Ufa, D; q, 1) is the free-particle 
propagator. In the path integral formulation of quantum mechanics, the propagator 
involves a sum over all possible space-time paths connecting the points (q, 1) and (x, D). 
The path integral expression for the free-particle propagator is a standard result, see 
for example. In the case at hand the free-particle propagator is given by 

12 



u 



1 



exp 



IX 



qi 



[2mv{D - l)] 3 / 2 1 [2u(D - 1) 
and so the wavefunction (IA.20 becomes 



[2mu(D - l)] 3 / 2 
Inserting the Madelung form 





r • i 1 2 i 
ix — q 


1 exp 


_2v{D - 1)_ 



^i(q)rf 3 q. 



1>i 



5i) 1/2 



exp 



-K 



of the initial wavefunction into flA.40 then yields 
1 



[2mu(D - 1)] 



3/2 



[l + <Ji(q)] 1/2 exp 



q| 



2(D-1) 



q 



rf 3 q. 



(A.3) 



(A.4) 



(A.5) 



(A.6) 



Recall that we are interested in the solution of the free-particle Schrodinger equation in 
the semi-classical limit v — > 0. The exponent appearing in (IA.6I) is complex and so the 
integrand will be a rapidly oscillating function of q in the limit u — > 0. The dominant 
contribution to the integral in flA.60 will then be from points where the phase varies 
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least rapidly with q, i.e. at stationary points. This method of evaluating the integral 
in flA.61) is referred to as the method of stationary phase. Before proceeding, we first 

J| by introducing a function S 



rewrite the exponent in 



S 



q| 



S(x, D; q, 1) of the form 
(A.7) 



2{D-1) 

Hereafter we will write S = S(q) for brevity since the integral in flA.61) is over q only. At 
stationary points q of the phase we then have S >m (q) = 0, where S >m denotes dS/dq m , 
leading to 

x = q+(D-l)s(q), (A.8) 

with the time-independent vector field s = s(q) defined by s = — V q 0i. In what follows 
we assume that, for a given value D of the linear growth factor, there is one and only 
one stationary point q satisfying ( I A. 81) for each x. This is tantamount to assuming the 
absence of shell-crossing. We now Taylor expand the function S about the stationary 
point q: 

-^3 3 

S(q) « S{q) + - ^ivm - q m )S, mn (q)(qn - q n ), (A.9) 

m=l n=l 

where S tTnn denotes d 2 S/dq m dq n and we have used the fact that 5, 
(1A.9I) into ( 1A.6I) and changing variables to p = q — q we obtain 

[i + <Mq)] 1/2 



(q) = 0. Inserting 



cxp 



■5(q) 



exp 



1 



-^p T M(q)p 



drp, 



(A.10) 



[2ttw{D - l)] 3 / 2 

where the superscript T denotes a transpose and M(q) is a complex symmetric matrix 
whose elements are M mn (q) = — iS )mn (q)/V. The three-dimensional Gaussian integral 
in flA.lOp can be evaluated using standard techniques: 



exp 



-p T M(q)p 



d 3 p 



(2vr) ; 
L^(q)J 



1/2 



(A.ll) 



where A^(q) is the determinant of M(q). In order for ( 1A.11I) to hold we require 
A^(q) 7^ 0; this condition is satisfied before shell-crossing occurs. To see this, first note 
that (IA.7D implies S >mn (q) = J mn (q)/(D - 1) where J mn (q) = S mn + (D - l)s m>n (q) 
and 5 mn is the Kronecker delta. The real symmetric matrix J(q) with elements J mn (q) 
is then simply the Jacobian matrix of the Zeldovich approximation evaluated at the 
stationary point q. It follows that 

i 



■M(q) 



(A.12) 



K^-i)] 3 ' 

where J{(\) is the determinant of the Jacobian matrix J(q). Since J{(\) is non-zero up 
until shell-crossing occurs, then so is A4(q). We can then substitute flATT2]) and flATll 
into OA-IPO to obtain 

"l + 5i(q) nl/2 



^(q) J 



cxp 



-S(q) 



(A. 13) 
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Comparing (1A.13j) with the Madelung transformation we see that, in the limit v 
the CDM density field 5 = <5(x, D) in the free-particle approximation is given by 



- m (A ' 14) 

and the velocity potential <p — <K X ; D) follows from (IA. T|) : 

= <Mq)-^OD-l)|s(q)| 2 , (A.15) 

where x and q satisfy (IA.80 . The expressions (1A.14D and (1A. 15[) are identical to those 
obtained in the Zeldovich approximation; see fflUj) and ffl3l) . Therefore we have shown 
that, prior to shell-crossing (i.e. while there is a one-to-one correspondence between x 
and q), the free-particle approximation reduces to the Zeldovich approximation in the 
semi-classical limit v — > 0. 



Appendix B. Perturbative solution of the one-dimensional Schrodinger 
equation with a time-dependent potential 

We have used the wave-mechanical approach to study the gravitational evolution of a 
plane-symmetric sinusoidal density perturbation in a static universe. The appropriate 
one-dimensional Schrodinger equation to be solved was 

dip ( v 2 d 2 
~dt ~ ~I& 2 



with an external gravitational potential $ = <&(x,t) of the form 



$ = D 



^-\{D-l)v 2 



(B.2) 



where = $i(q r ) is the initial gravitational potential, V[ = Vi(q) is the initial velocity 
field and Eulerian coordinates x are related to Lagrangian coordinates q by the Zeldovich 
mapping x = q+r(D — l)v\ with D = exp [(t — U)/t\. For an initial density perturbation 
<5i — $i(<l) °f the form 5; = <5 a cos (2-nq/d) we have 

^ = -(7^-X^ (B.3) 



(Si ~ 6?) ■ (B.4) 



2tit 

and 



Z7TT 



We now describe how to construct a second-order perturbative solution of the 
Schrodinger equation (IB. 1|) using TDPT. To begin with, note that (IB. II) is written in 
the position representation. It can be written in a more general form without referring 
to a particular basis as 

^ \m) = H(t) \m) , (b.5) 
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where \ip{t)) is a general ket representing the state of the physical system at time t 
and H{t) is the time-dependent total Hamiltonian. In order to solve (IB.5j) using PT 
it is convenient to split the total Hamiltonian according to H(t) = Hq + $(£) where 
H — -P 2 /2 is the free-particle Hamiltonian (with P the momentum operator) and $(£) 
is the time-dependent potential (IB. 21) in operator form. Before continuing to develop a 
perturbative solution of (IB. 51) we first describe how to solve the free-particle Schrddinger 
equation 

\m) = Ho \m) ■ (b.6) 

Appendix B.l. The free-particle Schrddinger equation 

The solution \ip(t)) of the general free-particle Schrddinger equation (IB. 61) is obtained 
from the initial state ket = \ip(ti)) via 

m)) = Uo(t,ti)\1>i), (B.7) 

where U (t,t{) is a unitary operator, known as the free-particle time-evolution operator, 
obeying 

uv-^o(Mi) =H U (tA)- (B.8) 
at 

The free-particle Hamiltonian H is independent of time and so we can directly integrate 
(IB~8|) to find that 



U (t,ti) = exp 



(B.9) 



The initial state ket can be expanded in terms of orthonormal eigenkets \n^) of 
H as 

|^i) = ^«n|n (0) ), (B.10) 

n 

where a n = (n^\ipi) and the eigenkets \n^} satisfy the general time-independent free- 
particle Schrddinger equation 

H \n®) = \n®) . (B.ll) 

In order to determine the eigenvalues En and the eigenkets |n^) we must solve the 
eigenvalue problem (IB. 111) . To do this it is convenient to first rewrite it in the position 
representation: 

where <pn^ = 4>n\x) are the eigenfunctions, defined by (j)^ = (x\n^). Introducing 
k\ = 2En^ /u 2 we find that the solutions to (IB.12I) are of the form 

= TJ/2 ex P( ifc nx), (B.13) 
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where k n = 2nn/L, n an integer, since we are considering a cubic volume of side length L 
equipped with periodic boundary conditions at each face. The pre-factor 1/L 3 / 2 comes 
from the requirement that the eigenf unctions are normalized. 

The solution of the free-particle Schrddinger equation (1B.6j) follows upon inserting 
the time-evolution operator (1B.9I) and the initial ket expansion (IB. 101) into (IB. 71) . In the 
position representation we find 



v> = £ 



a n exp 



-i{t - t{)E, 



(0) 



V 



do) 

'n i 



(B.14) 



(o) 



where if) = (x\if)(t)), E., 
expansion coefficients a r 



i >2 k 2 l /2 and the eigenfunctions 
(n^l^i) are found from 



do) 



are given by ( IB. 131) . The 



L 2 I <f> { n\qmq)dq, (B.15) 



o 



where ipi = (q\ipi) and an over-line denotes complex conjugation. As discussed 
previously, we divide the cubic volume into cells of side length d by setting L = Nd, 
N > an integer. We can then write ( IB. 151) as 



JV-l 



r(j+l)d 

i n = L 2 J2 / <Pn\q)^(q)dq. (B.16) 

j=0 J 3 d 



Performing a change of variable q' = q — jd leads to 

r d N-l 

a n = L 2 0i o) (g')^i(g / )^ , ^exp(-ijA; ri rf), (B.17) 
Jo j=0 

where we have used (IB. 13[) along with the fact that the initial Madelung wavefunction 
= (1 + 5 i ) 1 / 2 exp (—iipi/v) is periodic with period d. However, 

N-l 

^ exp {-\jk n d) = NS ntpN , (B.18) 

3=0 

where p is an integer, and thus the coefficients a n are given by 

d 



a n = NL 2 / $\q')^{q')dq' (B.19) 



o 



where the eigenfunctions (f>n are given by (IB.13j) . but with k n = 2nn/d. It follows from 
(E33D that the expansion coefficients are simply found by taking the Fourier transform 
of the initial wavefunction. 
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Appendix B.2. The Schrddinger equation with a time- dependent potential 

Our aim now is to perturbatively solve the Schrddinger equation ( lB.5fl with a time- 
dependent Hamiltonian of the form H(t) = Hq + <E>(t) where the potential is given by 
( 1B.2|) in the position representation. As before, the solution \ip(t)) of the Schrddinger 
equation is related to the initial ket via 

|</>(t)) = t>(Mi)|^i>, (B.20) 
where the time-evolution operator obeys 

iujtl{t^=H{t)U{t^. (B.21) 

The previously determined eigenkets \n^) of the free-particle Hamiltonian Ho form an 
orthonormal set and so we may expand the initial state ket as in flB.101) with the 
expansion coefficients a n = (n^°'\ipi) given by ( IB. 191) . In order to determine the state 



ket \ip(t)) at a time t, it is clear from (IB. 201) that we require an expression for the 
time-evolution operator U(t,tj). Since the Hamiltonian now depends explicitly on time, 
we can no longer simply integrate (IB. 211) to obtain an expression for the time-evolution 
operator as we did in the free-particle case. The strategy in such a situation is to seek 
an approximate expression for the time-evolution operator. We now describe how this 
may be achieved by using TDPT. 



Appendix B.2.1. Time- dependent perturbation theory In this work we have exclusively 
used a description of quantum dynamics known as the Schrddinger picture. However, 
TDPT solutions of the Schrddinger equation with a time-dependent potential are most 
easily constructed by appealing to another equivalent description of quantum dynamics 
known as the interaction picture. A general state ket in the interaction picture \ij){t)y' 
is related to the state ket in the Schrddinger picture \ip(t)) by 

mt)f ] = m,ti)mt)), (B.22) 

where the free-particle time-evolution operator U (t,t{) is given by ( IB. 91) and a dagger 
denotes the adjoint. The superscript (I) will be used throughout to denote quantities in 
the interaction picture. Note that = \ipi) and so the interaction and Schrddinger 

picture state kets initially coincide. In the interaction picture it is straightforward to 
show that a general state ket evolves according to: 

w±W))® = * m (t)m)) m , (B.23) 

where $®(t) = ?7q (t, ti)&(t)U (t, t ; ) is the perturbing potential in the interaction picture. 
It is apparent from (1B.23j) that the time-evolution of a state ket in the interaction picture 
is determined solely by ^^(t). 

The objective is to determine a perturbation expansion for the time-evolution 
operator U (t, t{) in the Schrddinger picture. This is best achieved by first finding a 
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perturbation expansion for the time-evolution operator in the interaction picture. In 
this picture a time-evolution operator U®(t,ti) can be defined by 

m) f) = um {tA)m m. (B.24) 

Substituting flR24p into (lB~23l) it is clear that U®(t,ti) obeys 

^4f> (I) (^ii) = ^ (I) Wf> (I) (t^i), (B.25) 
at 

where this equation must be solved subject to the initial condition U®(t{,ti) = I with 
/ the identity operator. Observe that equation flB.250 . together with the appropriate 
initial condition, is equivalent to the integral equation 

U® (t, ^ = 1-- f {t')U {1) (t f , U)dt'. (B.26) 
v Jti 

The integral equation (IB. 261) provides a convenient means of determining a perturbation 
expansion for U^(t,t{). By iteration we find that 

i>t i ft ft 1 

U^(t,ti) = I- - / ^'<f> (I) (t')--7 / dt' dt"&V(t')¥ l \t"), (B.27) 

v Jti v Jti Jti 

to second-order. To find the corresponding time-evolution operator U(t, t{) in the 
Schrddinger picture, first note that \ip(t))® = U^(t,ti)\ipi)^ = U^(t,ti)\ipi) since 
the state kets in the interaction and Schrddinger pictures coincide at t = U. Using the 
definition ( IB. 221) it follows that 

m)) = u (t,t i )uV(t,t i )m. (b.28) 

Upon comparison with (IB. 201) it is immediately clear that the time-evolution operator 
in the Schrddinger picture U(t,t[) is related to U®(t,ti) via 

f>(t,ti) = C>o(t,ti)t> (I) (t,ti). (B.29) 

Multiplying (lB~27l) by U (t,ti), inserting $W(t) = Ui(t,ti)$(t)U (t,ti) and using the 
following property of the time-evolution operator: Uo(t, ti)Uo(t' ', t{) = Uo{t, ti)Uo(ti, t') = 
Uo(t, if) then gives 

2 

U(t,t i ) = J2u j (t,t i ), (B.30) 

3=0 

to second order, where Uo(t,ti) is given by (1B.9|) . 

Ui{t,t i ) = -- [ dt'U {t,t')${t')U Q {t',t{), (B.31) 
v Jti 

and 

U 2 (t,ti) = -^f t dt 'f t c/t"f/o(t,* , )^(*')^o(t / ,0^(*")^o(t",ti)- (B.32) 
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The second-order expression for the time-evolution operator (1B.30|) can be inserted into 
( IB. 201) along with the initial ket expansion (IB. 101) to give the following approximate 
solution of the Schrddinger equation ( IB. 51) : 



where 



Mt)> = £W w (*)> 

3=0 



\^\t)) = J2*nU J (t,t i ) \n {0) ) 



(B.33) 



(B.34) 



with U (t,ti), C/i (t, *i) and U 2 (t 1 t i ) given by (lB~9l) . fjBT3l|) and (IB~32l) . respectively. 
In the position representation the zeroth-order term ip^ = (x\ip^(t)) is simply the 
solution (IB. 141) of the free-particle Schrddinger equation discussed earlier. Since the 
eigenfunctions 0^ of the free-particle Hamiltonian are independent of time, we may 
write the first-order term ip^ = (x\i/j^ (t)) as 



(1) 



(B.35) 



n m 



where 2 m ,n = 1 m ,n if) is defined by 



t 

1 WJ , = I dt'exp 



-i(t - t')E, 



(o) 

171 



V 



$ m , n (t')exp 



-i(t' - ti)K 



(o) 



V 



, (B.36) 



and the matrix elements $ m>n (t) = (m ( - - ) |$(t)|n^ - ) ) are given by 



tf(x)$(x,t')^ 0) (^)^- 



(B.37) 



In order to calculate the matrix elements $ m>n , we first use the Zeldovich approximation 
x = q + t[D — A)v\ to change to Lagrangian coordinates q. Inserting $ from ( 1B.2I) and 
using (1B.13|) then leads to 



L 



D 



H<l)-\{D-l)v!{q) 



exp (—ik p q) exp [—ik p r(D — l)uj(g)] 



xj(q,t')dq, 



(B.38) 



where k p = k m — k n and the Jacobian determinant J7" = dx/dq = 1 + t(D — l)dtij/dg. 
Upon substituting $i and v- t from (1B.3[) and ( IB. 41) we obtain, after a lengthy calculation, 



5„. ( d 



2 \ 2tyt 



D < a 



Jp{2pa) + - Jp+si^pa) 



s=±2 



cr 

T 



J p+S (2pa) 



(B.39) 



=±3 
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where p = m — n is an integer, a = a(t) is defined by a = S a (D — l)/2 and 

1 



Mx) 



2ty 



exp (iW) exp [— ix sin (0)]c?0 



(B.40) 



are Bessel functions of the first kind. In a similar manner we find that the second-order 
term tp^ = (x\ip( 2 \t)) is given by 



(2) 



■3 Y an Y Y fci,™,^' 



(B.41) 



m i 



where £j, m , n = K,^ m>n {t) is defined by 



,„.„ , r//' / dt" exp 
ti Jti 



(0) 



X$m,n(*") exp 



-i(£ , 

v $z, m (t>xp 

"_i( t " -t^' 



-i(f - t")E, 



(0) 
m 



7/ 



(B.42) 



and the matrix elements $/, m (t) = (i^|$(t)|m^) and $ m , n (t) = (m (0) |<l(t)|n (0) ) are 
obtained from (1B.39[) . 
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